home *** CD-ROM | disk | FTP | other *** search
/ MacHack 2000 / MacHack 2000.toast / pc / The Hacks / 3D QuickTime Dynamics / kSan Sources / kSanSimBoxes.c < prev    next >
Encoding:
C/C++ Source or Header  |  2000-06-24  |  11.6 KB  |  450 lines
  1. #include <fp.h>
  2.  
  3. #include "kSanSimulationPrivate.h"
  4. #include "OHApplicationPublic.h"
  5.  
  6.  short   resizeBoxHandle  (kSanSimulation *sim, long x, long y, long z) ;
  7. short   resizeBoxHandle  (kSanSimulation *sim, long x, long y, long z) 
  8. {  // slots is the number of particles per box we allocate in memory 
  9.     short err = noErr;
  10.     long i;
  11.     box *boxes;
  12.     long oldNumberOfBoxes = sim->howManyBoxes;
  13.     sim->boxDim.a = x;
  14.     sim->boxDim.b = y;
  15. #ifndef __2DVectors__
  16.     sim->boxDim.c = z;
  17. #endif
  18.     sim->howManyBoxes = x * y * z;    
  19.     
  20.     if (sim->howManyBoxes > oldNumberOfBoxes)
  21.     {
  22.         err = LCHSetHandleSize(&sim->boxArrayLCH,(sizeof(box) * (sim->howManyBoxes +2L))); 
  23.         if (!err)
  24.         {
  25.             
  26.             boxes =  simLockBoxes(sim); 
  27.             for (i = oldNumberOfBoxes+2; i < sim->howManyBoxes+2; ++i)
  28.             {
  29.                 initializeList (&(boxes[i].boxPartList));
  30.             }
  31.             simReleaseBoxes(sim);
  32.         }
  33.         else
  34.         {
  35.             sim->howManyBoxes = oldNumberOfBoxes;
  36.         }
  37.     }
  38.     else
  39.     {
  40.         boxes =  simLockBoxes(sim); 
  41.         for (i = sim->howManyBoxes+2 ; i < oldNumberOfBoxes+2; ++i)
  42.         {
  43.             disposeList (&(boxes[i].boxPartList));
  44.         }
  45.         simReleaseBoxes(sim);    
  46.         err = LCHSetHandleSize(&sim->boxArrayLCH,(sizeof(box) * (sim->howManyBoxes +2L))); 
  47.     }
  48.     return (err);
  49. }
  50. short   initBoxes  (kSanSimulation *sim) 
  51. {  // slots is the number of particles per box we allocate in memory 
  52.     short err = noErr;
  53.     long maxx,maxy,maxz;
  54.     long x[3],y[3],z[3];
  55.     long i,j,k,counter;
  56.     long xcount, ycount, zcount;
  57.     long howManyNBs;
  58.     long thisBoxLabel;
  59.     long thatBoxLabel;
  60.     Handle tempPadding = nilPointer;
  61.     box * thisBox;
  62.     box * boxes;
  63.     float dummy;
  64.  
  65. #ifdef __2DVectors__
  66. #pragma unused (zcount)
  67. #endif
  68.  
  69.     dummy = sim->simCell.cellDim.a/(sim->outerCutoff + 1.0); // implicit type conversion
  70.     maxx = dummy; 
  71.     dummy = sim->simCell.cellDim.b/(sim->outerCutoff + 1.0); 
  72.     maxy = dummy;
  73.  
  74. #ifndef __2DVectors__
  75.     dummy = sim->simCell.cellDim.c/(sim->outerCutoff + 1.0);
  76.     maxz = dummy;
  77. #else 
  78.     maxz = 1;
  79. #endif
  80.  
  81.     if (!(maxx>0)) maxx = 1;
  82.     if (!(maxy>0)) maxy = 1;
  83.     if (!(maxz>0)) maxz = 1;
  84.     
  85.     tempPadding = NewHandle(100000L);
  86.  
  87.     if (!err) err = resizeBoxHandle(sim, maxx, maxy, maxz);
  88.  
  89.     while ((err) && (sim->howManyBoxes > 9)) 
  90.     {
  91.         if (maxx > 3) maxx = maxx - 1;
  92.         if (maxy > 3) maxy = maxy - 1;
  93. #ifndef __2DVectors__
  94.         if (maxz > 3) maxz = maxz - 1;
  95. #endif 
  96.         err = resizeBoxHandle(sim, maxx, maxy, maxz);
  97.     }
  98.     if (err) 
  99.     {
  100.         doAlert("\pMemory error in box allocation");
  101.         sim->flags.error = kSanErrBoxMemory;
  102.     }
  103.     else
  104.     {
  105.         boxes = simLockBoxes(sim);
  106.         // initialize the which contains hyperspace particles
  107.         thisBox = &boxes[sim->howManyBoxes + 1];
  108.         thisBox->label = sim->howManyBoxes + 1;
  109.         thisBox->Xindex = -1;
  110.         thisBox->Yindex = -1;
  111.         thisBox->Zindex = -1;
  112.         for (howManyNBs = 0;howManyNBs<NUMNEIGHBORINGBOXES;++howManyNBs)
  113.         { // initialize the hyperspace box (particles not in cell)
  114.             thisBox->neighboringBoxes[howManyNBs] = 0L;  // setting all the neighbor boxes to 0 means that 
  115.               // particles in box zero only interact with particles in box 0 :: i.e. no particles
  116.               // particles must be initialized to box (sim->howManyBoxes + 1);
  117.         }
  118.         // initialize the box which never contains anything
  119.         counter = 0;
  120.         thisBox =  boxes;
  121.         thisBox->label = counter;
  122.         thisBox->Xindex = -1;
  123.         thisBox->Yindex = -1;
  124.         thisBox->Zindex = -1;
  125.         for (howManyNBs = 0;howManyNBs<NUMNEIGHBORINGBOXES;++howManyNBs)
  126.         { // initialize the zero box (not in cell)
  127.             thisBox->neighboringBoxes[howManyNBs] = 0L;  // setting all the neighbor boxes to 0 means that 
  128.               // particles in box zero only interact with particles in box -1 :: i.e. no particles
  129.               // particles must be initialized to box 0;
  130.         }
  131.         ++ counter;
  132.         // label neighbors of all other boxes
  133.         for (i=0;i<maxx;++i)
  134.         {
  135.             if (i==0) x[0] = maxx-1;
  136.             else x[0] = i-1;
  137.             x[1] = i;
  138.             if (i==(maxx-1)) x[2] = 0;
  139.             else x[2] = i+1;
  140.             for (j=0;j<maxy;++j)
  141.             {
  142.                 if (j==0) y[0] = maxy-1;
  143.                 else y[0] = j-1;
  144.                 y[1] = j;
  145.                 if (j==(maxy-1)) y[2] = 0;
  146.                 else y[2] = j+1;
  147.                 for (k=0;k<maxz;++k)
  148.                 {
  149.                     thisBox = &boxes[counter];
  150.                     //thisBox =  (boxes + ((counter) * (sim->boxRecordSize)));
  151.                     thisBox->label = counter;
  152.                     thisBox->Xindex = i;
  153.                     thisBox->Yindex = j;
  154.                     thisBox->Zindex = k;
  155.  
  156.                     if (k==0) z[0] = maxz-1;
  157.                     else z[0] = k-1;
  158.                     z[1] = k;
  159.                     if (k==(maxz-1)) z[2] = 0;
  160.                     else z[2] = k+1;
  161.  
  162.                     thisBoxLabel = maxz*maxy*i + maxz*j + k + 1L;
  163.                     if (thisBoxLabel != counter) doAlert("\pBox counting error. Unknown cause. Exiting");
  164.                     howManyNBs = 0;
  165.                     for (xcount=0;xcount<3;++xcount)
  166.                     {
  167.                         for (ycount=0;ycount<3;++ycount)
  168.                         {
  169. #ifndef __2DVectors__
  170.                             for (zcount=0;zcount<3;++zcount)
  171.                             {
  172.                                 thatBoxLabel = maxz*maxy*x[xcount] + maxz*y[ycount] + z[zcount] + 1L; // include same box
  173.                                 thisBox->neighboringBoxes[howManyNBs] = thatBoxLabel;
  174.                                 ++howManyNBs;
  175.                             }
  176. #else
  177.                                 thatBoxLabel =  maxy*x[xcount] +  y[ycount]   + 1L; // include same box
  178.                                 thisBox->neighboringBoxes[howManyNBs] = thatBoxLabel;
  179.                                 ++howManyNBs;
  180. #endif
  181.                         }
  182.                     }    
  183.                     ++counter;
  184.                 }
  185.             }
  186.         }
  187.         for (thisBoxLabel=0; thisBoxLabel <= sim->howManyBoxes + 1L ; ++thisBoxLabel)
  188.         {
  189.             thisBox = &boxes[thisBoxLabel];
  190.             //thisBox =    (boxes + (thisBoxLabel * (sim->boxRecordSize)));
  191.             thisBox->boxPartList.howManyItems = 0;  // clear box contents
  192.         }
  193.         simReleaseBoxes(sim);
  194.         OHDisposeHandle(&tempPadding);
  195.         sim->counters.neighbor = sim->countThreshholds.neighbor + 1L;    // force a complete get neighbors
  196.     }
  197.  
  198.     return (err);
  199. }
  200. short  setBoxForPartList(kSanSimulation *sim, OHList *partList)
  201. {
  202.     short err  = noErr;
  203.     particle *unusedBasePart;
  204.     box * unusedBoxes;
  205.     // these calls should speed up the function by preventing toolbox calls to
  206.     // LockHandle and UnlockHandle inside sim->setBoxForPart()
  207.     unusedBasePart = simLockParts(sim);
  208.     unusedBoxes = simLockBoxes(sim); 
  209.     // end unused stuff
  210.     
  211.     err =  OHListEveryItemAsSecondArg( partList, sim->setBoxForPart, sim);
  212.     
  213.     simReleaseBoxes(sim);
  214.     simReleaseParts(sim);
  215.     return (err);
  216. }
  217.  
  218. long  determineBoxForPartTriclinic(kSanSimulation *sim,long partArrIndex)
  219. {
  220.     short err  = noErr;
  221.     long whichBox;
  222. //    particle * thisPart;
  223.     particle * basePart;
  224.     short whichX;
  225.     short whichY;
  226.     short whichZ;
  227.     double XCutOff;
  228.     double YCutOff;
  229.     double ZCutOff;
  230.     double inverseCellXDim;
  231.     double inverseCellYDim;
  232.     double inverseCellZDim;
  233.     kSanVector  adjPos;
  234.     
  235. #ifdef __2DVectors__
  236.    #pragma unused(whichZ)
  237.    #pragma unused(ZCutOff)
  238.    #pragma unused(inverseCellZDim)
  239. #endif
  240.     
  241.     basePart = simLockParts(sim);
  242. //    thisPart = &(basePart[partArrIndex]);
  243.  
  244.     XCutOff = sim->simCell.cellDim.a *0.5;
  245.     YCutOff = sim->simCell.cellDim.b *0.5;
  246. #ifndef __2DVectors__   
  247.     ZCutOff = sim->simCell.cellDim.c *0.5;
  248. #endif
  249.     inverseCellXDim = sim->boxDim.a / (sim->simCell.cellDim.a + 0.01); // these small number give a little leeway
  250.     inverseCellYDim = sim->boxDim.b / (sim->simCell.cellDim.b + 0.01);
  251. #ifndef __2DVectors__  
  252.     inverseCellZDim = sim->boxDim.c / (sim->simCell.cellDim.c + 0.01);
  253. #endif    
  254.     err = getAdjustedPositionWithWrap(&sim->simCell,  &basePart->basePosition[partArrIndex], &adjPos);    
  255.  
  256.     if (!err) 
  257.     {
  258.         whichX = trunc((adjPos.a + XCutOff) * inverseCellXDim);
  259.         whichY = trunc((adjPos.b + YCutOff) * inverseCellYDim);
  260. #ifndef __2DVectors__  
  261.         whichZ = trunc((adjPos.c + ZCutOff) * inverseCellZDim);
  262.         whichBox = sim->boxDim.c*sim->boxDim.b*whichX + sim->boxDim.c*whichY + whichZ + 1;
  263. #else
  264.         whichBox =  sim->boxDim.b*whichX +  whichY  + 1; 
  265. #endif
  266.         if (whichBox > sim->howManyBoxes) whichBox = -1;
  267.     }
  268.     else
  269.     {
  270.         whichBox = -1;
  271.     }
  272.     simReleaseParts(sim);
  273.     
  274.     return (whichBox);
  275. }
  276.  
  277. short  setBoxForPartTriclinic(kSanSimulation *sim,long partArrIndex)
  278. {
  279.     short err  = noErr;
  280.     box * boxes;
  281.     particle *basePart;
  282.     long whichBox;
  283.     long oldBox;
  284.     
  285.     basePart = simLockParts(sim);
  286.     boxes = simLockBoxes(sim);
  287.     whichBox =   determineBoxForPartTriclinic(sim, partArrIndex);
  288.  
  289.      oldBox = BPNthBoxIndex(basePart, partArrIndex);
  290.     if ((oldBox != whichBox) && (oldBox != -1))
  291.     {
  292.         removeItem(&boxes[oldBox].boxPartList, partArrIndex); 
  293.         err = addPartToBox(sim, basePart, partArrIndex, whichBox, boxes);
  294.         switch (err)
  295.         {
  296.             case noErr:
  297.             default:
  298.                 break;
  299.             case kSanErrUnallocatedBox:
  300.                 doAlert("\p unallocated neighbor box error. Please report this error and how you arrived here!");
  301.                 break;                
  302.         }
  303.     }
  304.     simReleaseParts(sim);
  305.     simReleaseBoxes(sim);
  306.     return (err);
  307. }
  308. long  determineBoxForPartOrthorhombic(kSanSimulation *sim,long partArrIndex)
  309. {
  310.     long err  = noErr;
  311.     long whichBox;
  312.     //particle * thisPart;
  313.     kSanVector  *thisPos;
  314.     particle * basePart;
  315.     short whichX;
  316.     short whichY;
  317.     short whichZ;
  318.     char nanFlag = false;
  319.     char outOfBox = false;
  320.     double XCutOff;
  321.     double YCutOff;
  322.     double ZCutOff;
  323.     double inverseCellXDim;
  324.     double inverseCellYDim;
  325.     double inverseCellZDim;
  326.  
  327. #ifdef __2DVectors__
  328.    #pragma unused(whichZ)
  329.    #pragma unused(ZCutOff)
  330.    #pragma unused(inverseCellZDim)
  331. #endif
  332.  
  333.     basePart = simLockParts(sim);
  334.     //thisPart = &(basePart[partArrIndex]);
  335.     thisPos = &basePart->basePosition[partArrIndex];
  336.     XCutOff = sim->simCell.cellDim.a *0.5;
  337.     YCutOff = sim->simCell.cellDim.b *0.5;
  338.  
  339. #ifndef __2DVectors__  
  340.     ZCutOff = sim->simCell.cellDim.c *0.5;
  341. #endif
  342.  
  343.     if ((thisPos->a <= XCutOff) && (thisPos->a >= ( -XCutOff))) {}
  344.     else 
  345.     { 
  346.         if (thisPos->a > XCutOff)
  347.         {
  348.             thisPos->a -= sim->simCell.cellDim.a;
  349.             if(thisPos->a > XCutOff) outOfBox = true;
  350.         }
  351.         else  if (thisPos->a < (-XCutOff))
  352.         {
  353.             thisPos->a += sim->simCell.cellDim.a;
  354.             if (thisPos->a < (-XCutOff)) outOfBox = true;
  355.         }
  356.         else  outOfBox = true;
  357.      }
  358.     
  359.     if ((thisPos->b <= YCutOff) && (thisPos->b >= ( -YCutOff))) {}
  360.     else 
  361.     { 
  362.         if (thisPos->b > YCutOff)
  363.         {
  364.             thisPos->b -= sim->simCell.cellDim.b;
  365.             if(thisPos->b > YCutOff) outOfBox = true;
  366.         }
  367.         else  if (thisPos->b < (-YCutOff)) 
  368.         {
  369.             thisPos->b += sim->simCell.cellDim.b;
  370.             if (thisPos->b < (-YCutOff)) outOfBox = true;
  371.         }
  372.         else  outOfBox = true;
  373.      }
  374.  
  375. #ifndef __2DVectors__  
  376.     if ((thisPos->c <= ZCutOff) && (thisPos->c >= ( -ZCutOff))) {}
  377.     else 
  378.     { 
  379.         if (thisPos->c > ZCutOff)
  380.         {
  381.             thisPos->c -= sim->simCell.cellDim.c;
  382.             if(thisPos->c > ZCutOff) outOfBox = true;
  383.         }
  384.         else   if (thisPos->c < (-ZCutOff))
  385.         {
  386.             thisPos->c += sim->simCell.cellDim.c;
  387.             if (thisPos->c < (-ZCutOff)) outOfBox = true;
  388.         }
  389.         else  outOfBox = true;
  390.      }
  391. #endif  
  392.      
  393.     if (!outOfBox)
  394.     {
  395.         inverseCellXDim = (sim->boxDim.a - 0.001) * sim->simCell.invCellDim.a  ; // these small number give a little leeway
  396.         inverseCellYDim = (sim->boxDim.b - 0.001) * sim->simCell.invCellDim.b  ; // these small number give a little leeway
  397. #ifndef __2DVectors__  
  398.         inverseCellZDim = (sim->boxDim.c - 0.001) * sim->simCell.invCellDim.c  ; // these small number give a little leeway
  399. #endif
  400.  
  401.         whichX = trunc((thisPos->a + XCutOff) * inverseCellXDim);
  402.         whichY = trunc((thisPos->b + YCutOff) * inverseCellYDim);
  403. #ifndef __2DVectors__  
  404.         whichZ = trunc((thisPos->c + ZCutOff) * inverseCellZDim);
  405.         whichBox = sim->boxDim.c*sim->boxDim.b*whichX + sim->boxDim.c*whichY + whichZ + 1;
  406. #else
  407.         whichBox =  sim->boxDim.b*whichX +  whichY   + 1;
  408. #endif
  409.         if ((whichBox < 0) || (whichBox > sim->howManyBoxes)) whichBox = -1;
  410.     }
  411.     else whichBox = -1;
  412.     
  413.     simReleaseParts(sim);
  414.     return (whichBox);
  415. }
  416.  
  417. short  setBoxForPartOrthorhombic(kSanSimulation *sim,long partArrIndex)
  418. {
  419.     long err  = noErr;
  420.     box * boxes;
  421.     long oldBox;
  422.     long whichBox;
  423.     particle *basePart;
  424.     
  425.  
  426.     basePart = simLockParts(sim);
  427.     boxes = simLockBoxes(sim);
  428.     
  429.     whichBox =  determineBoxForPartOrthorhombic(sim, partArrIndex);
  430.  
  431.      oldBox = BPNthBoxIndex(basePart, partArrIndex);
  432.     if ((oldBox != whichBox) && (oldBox != -1))
  433.     {
  434.         removeItem(&boxes[oldBox].boxPartList, partArrIndex); 
  435.         err = addPartToBox(sim, basePart, partArrIndex, whichBox, boxes);
  436.         switch (err)
  437.         {
  438.             case noErr:
  439.                 break;
  440.             default:
  441.             case kSanErrUnallocatedBox:
  442.                 OHSystemBeep();
  443.                 break;                
  444.         }
  445.     }
  446.     simReleaseParts(sim);
  447.     simReleaseBoxes(sim);
  448.     return (err);
  449. }
  450.